CONTINUOUS GALERKIN FINITE ELEMENT METHODS FOR 
HYPERBOLIC INTEGRO-DIFFERENTIAL EQUATIONS 



FARDIN SAEDPANAH 



Abstract. A hyperbolic integro-differontial equation is considered, as a model 
problem, where the convolution kernel is assumed to be either smooth or no 
worse than weakly singular. Well-posedness of the problem is studied in the 
context of semigroup of linear operators, and regularity of any order is proved 
for smooth kernels. Energy method is used to prove optimal order a priori error 
estimates for the finite element spatial semidiscrete problem. A continuous 
space-time finite element method of order one is formulated for the problem. 
Stability of the discrete dual problem is proved, that is used to obtain optimal 
order a priori estimates via duality arguments. The theory is illustrated by an 
example. 



1. Introduction 

We consider, for any fixed T > 0, a hyperbolic type integro-differential equation 
of the form 

(1.1) il + Au- j fC(t-s)Au(s)ds = f, te(0,T), with u(0) = u°, w(0) = u\ 

Jo 

(we use '•' to denote 'J^') where A is a self-adjoint, positive definite, uniformly 
elliptic second order operator on a Hilbert space. The kernel 1C is considered to be 
either smooth (exponential), or no worse than weakly singular, and in both cases 
with the properties that 

(1.2) /C>0, £(t)<0, \\JC\\ Ll (n+) = « < 1. 

This kind of problems arise e.g., in the thoery of linear and fractional order 
viscoelasticity. For examples and applications of this type of problems see, e.g., 
[13], [7], and references therein. 

For our analysis, we define a function £ by 

f-t />oo 

(1.3) £(t) = /s- / JC(s)ds= / JC(s)ds, 



and, having (1.2 1, it is easy to see that 

(1.4) A£(*) = -£(*) < 0, $(0) = «, Hm£(t)=0, <£(*)<«. 

Hence, £ is a completely monotone function, since 

(-l)'Dfc(t) > 0, te (0,00), j =0,1, 2, 
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and consequently £ € ii,/ O c[0, oo) is a positive type kernel, that is, for any T > 
and 0eC([O,T]), 

(1.5) / [ £(t-s)(f>(t)(f>(s)dsdt>0. 

Jo Jo 

From the extensive literature on theoritical and numerical analysis for partial 
differential equations with memory, we mention [T3], [7J, [2], [TO], [2], and their 
references. 

The fractional order kernels, such as Mittag-Leffler type kernels in fractional vis- 
coelasticity, interpolate between smooth (exponential) kernels and weakly singular 
kernels, that are singular at origin but integrable on finite time intervals (0,T), 
for any T > 0, see [TS] and references therein. This is the reason for considering 
problem (1.1 1 with convolution kernels satisfying (1.2). 

In [7] well-posedness of a problem, similar to (1.1) with a Mittag-Leffler type 
kernel, was studied in the framework of the linear semigroup theory. Here we first 
extend the theory to prove higher regularity of the solution for more smooth kernels, 
such that a priori error estimates are fulfilled. We prove L 00 (L2) optimal order a 
priori error estimate, by energy methods, for finite clement spatial semidiscrete 
approximate solution. This provides an alternative proof to what we presented in 
[7j, and is straightforward. The continuous space-time finite element method of 
order one, cG(l)cG(l), is used to formulate the fully dicrete problem. A similar 
method has been applied to the wave equation in [5], where adaptive methods 
based on dual weighted residual (DWR) method has been studied. An energy 
identity is proved for the discrete dual problem, using the positive type auxiliary 
function £. This is then used to prove L^L?) and L^H 1 ) optimal order a priori 
error estimates by duality. This and [14] , where a posteriori error analysis of this 
method has been studied via duality, complete the error analysis of this method for 
model problems similar to (1.1). 

The present work also extend previous works, e.g., [2J, [T], [TS], on quasi-static 
fractional order viscoelsticity (u rs 0) to the dynamic case. Spatial finite element 
approximation of integro-differcntial equations similar to (1.1) have been studied 
in [3] and [5], however, for optimal order L 00 {L2) a priori error estimate for the 
solution u, they require one extra time derivative regularity of the solution. A 
dynamic model for viscoelasticity based on internal variables is studied in [T3 . The 
memory term generates a growing amount of data that has to be stored and used in 
each time step. This can be dealt with by introducing "sparse quadrature" in the 
convolution term |19j . For a different approach based on "convolution quadrature" , 
see [T7]. However, we should note that this is not an issue for exponentially decaying 
memory kernels, in linear viscoelasticity, that are represented as a Prony series. In 
this case recurrence relationships can be derived which means recurrence formula 
are used for history updating, see [18] and [9] for more details. In practice, the 
global regularity needed for a priori error analysis is not present, e.g., due to the 
mixed boundary conditions, that calls for adaptive methods based on a posteriori 
error analysis. We plan to address these issues in future work. 

In the sequel, in §2, well-posedness of the problem is proved and high regularity 
of the solution of the problem with smooth kernels is verified. In §3, the spatial 
finite element discretization is studied and, using energy method, optimal order a 
priori error estimates are proved. The continuous space-time finite element method 
of order one is applied to the problem in §4, and stability estimates for the discrete 
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dual problem are obtained. These are then used to prove optimal order a priori 
error estimates in §5 by duality. Finally, in §6, we illustrate the theory by a simple 
example. 



2. WELL-POSEDNESS AND REGULARITY 

We use the semigroup theory of linear operators to show that there is a unique 



solution of (1.1 1, and we prove that under appropriate assumptions on the data 
we get higher regularity of the solution. In §2.1 we quote the main framework 
from [7], to prove existence and uniqueness, to be complete. Here we restrict to 
pure homogeneous Dirichlet boundary condition, though the presented framework 
applies also to mixed homogeneous Dirichlet-Neumann boundary conditions. But it 
does not admit mixed homogeneuos Dirichlet nonhomogencous Neumann boundary 
conditions, and this case has been studied in [15] for a more general problem, by 
means of Galerkin approximation method. Then in §2.2 we extend the semigroup 
framework to prove regularity of any order for models with smooth kernels. To this 
end, we specialize to the homogeneous Dirichlet boundary condition. 

2.1. Existence and uniqueness. We let f2 C R d , be a bounded convex domain 
with smooth boundary dQ,. In order to describe the spatial regularity of functions, 
we recall the usual Sobolev spaces H s = H s (Q) d with the corresponding norms 
and inner products, and we denote H = H° = L 2 (Q) d , V = H^(fl) d . We equip 
V with the energy inner product a(u,v) = (Au,v) and norm \\v\\y = a(v,v). We 
recall that A is a selfadjoint, positive definite, unbounded linear operator, with 
T>(A) = H 2 n V, and we use the norms \\v\\ s — \\A s / 2 v\\. We note that with mixed 
homogeneous Dirichlet-Neumann boundary conditions, we have 

V = {v e H 1 : v = on Dirichlet boundary}. 

We extend u by u(t) = h(t) for t < with h to be chosen. By adding — J_ K(t— 
s)Ah(s) ds to both sides of ( |1.1[ ), changing the variables in the convolution terms 
and defining w(t, s) = u{t) — u(t — s), we get 

/•oo />oo 

(2.1) u(t) + (l- K )Au{t) + J fC(s)Aw(t,s)ds = /(*) - J tC(s)Ah(t - s) ds, 
where, we recall that ||/C||^ 1 (r+) = n < 1. For latter use, we note that equation 



(1.1 ) can be retained from (2.1 1 by backward calculations. 

For a given integer number r > 0, we use the Taylor expansion of order r of the 
solution u at t = to define the extension u(t) — h r (t) for t < 0. That is, we set 

T 4-71 

(2.2) u(t) = h r (t) = V -ru n (0), t < 0, 

n=0 

where we use the notation u n (t) = u n (t, ■) = j^u(t, •), with u°(t) = u(t). 



Now we reformulate the model problem (1.1) to an abstract Cauchy problem 



First, we choose r = in (2.2), that is ho(t) — u°, and for the initial data we 
assume that u° € T3{A) and u l e V. Therefore, from (2.1), we have 



(2.3) u{t) + (I - k)Au(£) + / K,{s)Aw{t,s)ds = f(t) - Au° K.(s) ds, 

Jo Jt 
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where, 



w(t, s) 



u(t)-u(t-s), se[0,t], 
u(t) — u , s G [t, oo). 



Then we write (2.3), together with the initial conditions, as an abstract Cauchy 
problem and prove well-posedness. 

We set v = ii and define the Hilbert spaces 

W = L 2iK (R+; V) = [w : \\w\\ 2 w = J K{s) \\w(s) \\ 2 V ds < co}, 

Z = V x H x W = jz = (u,v, w) : \\zf z = (1 - k)||m||v + IMP + Hlw < °°}- 
We also define the linear operator A on Z such that, for z = (it, v, w), 

- Dw 



Az = (v, -A(jl - k)u + J K,{s)w(s)ds S j 



(2.4) 



with domain of definition 
V{A) = \{u,v,w) e Z : v e V, (1 - K ) u + IC{s)w(s)ds e V(A),w G P(-D)|. 

Here = £w with P(D) = {w £ W : Dw € W and w(0) = 0}. 

Therefore, a solution of (1.1) satisfies the system of delay differential equations, 
for t € (0,T), 

ii(t) = v, 

v(t) = -A((l - K)u(t) + J K{a)w(t,a)daj+f{t)-Au° J tC(s) ds, 

w(t, s) = v(t) — Dw(t, s), s G (0, oo). 
This can be writen as the abstract Cauchy problem 

z(t) = Az(t) + F(t), t€(p,T), 
z(0) = z°, 

where F(t) = (0, f(t) - Au° J™K,{s) ds,0) and z° = (u ^ 1 ^), since 
w{0, s) = u(0) - u(-s) = u(0) - h(-s) = u° - u° = 0. 

We note that w(t,0) = u(t) - u{t) = 0, so that w(t, ■) € V(D). 

We quote from [3 Theorem 2.2], that A generates a Co-semigroup of cotractions 
on Z. 

Corollary 1. The linear operator A is an infinitesimal generator of a Co-semigroup 
e t_4 of contractions on the Hilbert space Z . 

Now, we look for a strong solution of the initial value problem ( |2.4[ ), that is, a 
function z which is differentiable a.e. on [0,T] with z € Li((0,T); Z), if z(0) = z°, 
z(t) G T>(A), and z(t) = Az(t) + F(t) a.e. on [0,T]. 

Recalling the assumptions u° G F>{A) and m 1 G V, we know that if z = (u, u, to) 
be a strong solution of the abstract Cauchy problem (2.4) with z° — (vP^u 1 ^), 
then u is a solution of (1.1) by [7, Lemma 2.1]. Hence, to prove that there is a 
unique solution for (1.1), we need to prove that there is a unique strong solution 
for ([2T|. This has been proved in Theorem 2.2], if / : [0,T] -> H is Lipschitz 
continuous, using the fact that the linear operator A generates a Co-semigroup 
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of contractions on Z. Moreover, for some C — C(k,T), we have the regularity 
estimate, for t £ [0, T], 



\\u(t)\\ v + \\u(t)\\ < c(\\Au°\\ + Hu 1 !! + f ll/H da) 



(2.5) 



2.2. High order regularity. In order to prove higher regularity of order r + 1 
( r > 1), we assume that the bounded domain fl is convex, and we specialize to the 
homogeneous Dirichlet boundary condition. Hence, the elliptic regularity estimate 
holds, that is 



(2.6) 



|u|| 2 < C\\Au\\, u£V(A) = H 2 nV. 



We note that the case r — is the choice for (2.5 1. We substitute h r (t) from 



(2.2), with r > 1, in (2.1 ). Then, differentiating Jp and using the notation u r {t) — 
§£Fu(t), we have 

il r (t) + {l- K )Au r (t)+ I IC(s)Aw r (t,s)ds 
= f(t)-A 9 

(2.7) 



o 



dr 

r-l 



poo r 

/ WE- 

Jt n=0 



f r (t) +AY, u n (0)JC r - n - 1 {t) - Au r {0)£,(t) 



n=0 



= ■ fit), 

with the initial data u r (0), u r+1 (0). 

Recalling the initial data u(0) = u° and m 1 (0) = it 1 , from (1.1 1, we have w 2 (0) = 
/(0) — Au°. To obtain u" l (0), to > 3, we differentiate g t „.~2 of equation ( |1.1[ ), and 
we have 



(2.; 



+ Au m - 2 (t) - I K m -\t- s)Au(s)ds 

m— 3 



/ m ~ 2 (t ) + A u " (t)/C m -"- 3 (0) , t e (0, T), 



n=0 



that, with t = 0, implies the initial condition 

m— 3 



(2.9) 



l (0) = / m ~ 2 (0) - Au m - 2 (0) + Au n (0)lC 



m—n—S 



(0), to>3. 



Throughout, obviously any sum Xm=j ^ s supposed to be suppressed from the 
formulas, when i > j. 

Remark 1. We note that, if we assume K £ W™~ 2 (0, T), then K £ C m " 3 [0, T] by 
Sobolev inequality, and therefore u m (0) is well-defined. 



Remark 2. One can show, by induction and the fact that by (2.6) 
(2.10) HI < IWIh* < C\\Av\\, foxveV(A), 
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we have [m = 0, 1, k = 1, 2, • • • ), 



|A m u 2 *(0)| < c(|A 



1^' 



m+fc-l u l| 



(2.11) 



(2.12) 



fe-l fc-2 

^ | A m+ij2fc-2j-2( )| + ^ | A m+ J/ 2fe-2 J -3( )| 



|yl m it 2fc+1 (0)| < c(\A m+k u°\ + |A m +V| 



fc-1 fe-l 
+ ^ | yl » l +J / 2fc-2 3 -2 (0 )| + ^ \A m+j f k - 2j - 1 (0)\ 

3=1 3=0 



Now we note that, in (2.7), we have 

' u r {t)~u r {t- a), «e[0,t], 



iu r (i,s) 



li r (t) -U r (0), 



s e [t, oo), 



so that w r (t,0) = 0. Therefore, considering continuty of w r , we have w r £ T>(D). 

Then, in the same way as in the previous section, with v r — u r , we can reformu- 
late ( |2.7[ ), with z r — (u r ,v r ,w r ), as the abstract Cauchy problem 

z r (t) = Az r (t)+F r (t), 0<t<T, 

(2.13) 

V ' z r (0) = z r '°, 

where F r (t) = (0,/ r (i),0) and z r '° = (u r (0), w r+1 (0), 0), since w r (0,s) = u r {0) - 
u r (0) = 0. 

In particular, for r = 1, we have 



F 



\t) = (o > / 1 (t) + A«°K:(i)-A«^(t), 



with initial data z 1 - = (u 1 , u 2 (0), 0) = (u 1 , /(0) - 0). 

Now, we need to show that from a strong solution of the abstract Cauchy problem 
(2.13), for r > 1, we get a solution of the main problem (1.1). Therfore we should 
prove that the abstract Cauchy problem (2.13) has a unique strong solution, under 
certain conditions on the data. The proof is by induction, and therefore we recall 
some facts from [7], for r = 1. 

Lemma 1. Let z 1 — (u ,v ,W ) be a strong solution of the abstract Cauchy problem 
( |2.13| with z lfi = (n\ it 2 (0),0). Thenu(t) = u° + J^u 1 (s)ds is a solution of JlTJ. 



Theorem 1. There is a unique solution u = u(t) of (1.1) if K € W^M" 1 ") with 
I|£|Ivk 1 1 (R+) < !; /(0) - -4w° G V, m 1 G TJ(A), and f : [Q,T\ -> H is Lipschitz 
continuous. Moreover, for some C = C(k,T), we have the regularity estimate, for 

te[o,T], 

ft 



(2.14) 



IK*)||v + ||«(*)|| < C (\\Au°\\ + HAu 1 !! + ||/(0)|| 



\f\\ds 



Proof. There exists a unique strong solution z 1 — (it 1 , v 1 , w 1 ) for ( |2.13[ ), with r = 1, 
by Theorem 2.4]. Hence, the proof is complete by Lemma [l] □ 



Lemma 2. Let z r 



(u r , v r , w r ), for r > 1, be a strong solution of the abstract 
Cauchy problem pT3| with z r <° = (w r (0), u r+1 (0), 0). Thenu{t) = u° + J*u 1 (s)ds 
is a solution of (1.1). 
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Proof. The proof is by induction. The case r = 1 follows from Theorem [T] 

Now, we assume that the lemma valids for some r > 2, and we prove that it 
holds also for r + 1. To this end, we show that if z r+1 — (u r+1 ,v r+1 ,w r+1 ) be 
a strong solution of ( |2.13[ ) (for r + 1) with z r+1 '° = (it r+1 (0), u r+2 {0), 0), then 
z r = (u r ,v r ,w r ) is a strong solution of ( |2.13| with z r >° = (w r (0), it r+1 (0), 0), that 
completes the proof by induction assumption. 

Since z r+1 (t) = Az r+1 (t) + F r+1 {t) a.e. on [0,T], we have, for t e (0,T), 



v r+i (t) = -A[(l - K )u r+L (t) + J K(s)w r+1 (t, s) ds) + 

w r+1 (t,s) =v r+1 (t)- Dw r+1 (t,s), s£ (0,oo). 

The first and the third equation implies that w r+1 satisfies the first order partial 
differential equation 

w r t +1 + w r s +1 =u t 

This, with w r+1 (t 7 0) = 0, w r+1 (0,s) = 0, has the unique solution w r+1 (t, s) 
u r+1 (t) — u r+1 (t — s), that implies, by integration with respect to t, 

w r (t,s)=u r (t)-u r (t-s)= f w r+1 (T,s)dT. 



,r+l 
H 



From the first and the second equations we obtain equation (2.7) with r + 1, that 
is obtained from equation (2.1) by differentiating J^. We recall that equations 
( 1.1 ) and (2.1 ) are equivalent, that implies equivalence of equations (2.7 1 and (2. 



Therefore u also satisfies (2.8) with r + 1. Then, integrating with respect to t, we 
have, for t G (0,T), 



u r {t) - u r {0) + Au r {t) - Au r (0) - 



(2.15) 



rC r+i (r - s)Au(s) ds dr 



"'O 

t 



f{t) ~f r (0)+Aj2 [ u n (r) drK r ~ n 

n=0 Jo 



(0). 



Now, we need to show that (2.15) implies ( |2.8[ ). We note that 



"'O 



K. r+L (T-s)Au(s) ds dr 



AC r+i (r- s)Au(s) dr ds 



Js 
t 



IC r (t- s)Au{s) ds- A/C7(0) / u{s) ds, 



and 



A V f u n {r) drlC r - n (0) 

= a( Cu{t) dr/C r (0) + V (« n_1 (*) - " n_1 (0))/C r -"(0)) 
V^ „=i / 

„t r—1 r—1 

= A/C r (0) / u(t) dT + AY,u n (t)K. r - n - 1 (0)-AY,u n (0)K. r - n - 1 {0). 

•'° n=0 n=0 

Using these and (2.9) in ( 2.15| ) we conclude (2.8), that is equivalent to (2.7). 
This means that, z r — (u r ,v r ,w r ) is a strong solution of (2.13) with z r ' u = 
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(u r (0), ii r+1 (0), 0). Hence, by induction assumption, u(t) = u° + J^u 1 (s)ds is 
a solution of (1.1 1, and this completes the proof. □ 

In the next theorem we find the circumstances under which there is a unique 
strong solution of the abstract Cauchy problem (2.13), that by Lemma [2] implies 
existence of a unique solution of ( |1.1[ ) with higher regularity. We also obtain regu- 
larity estimates, which are extensions of (2.5) and ( |2.14 ). 

We note that, recalling Remark [I] and having the assumptions from the next 
theorem, the calculations in the proof of Lemma [2] make sense. 

Theorem 2. For a given integer number r > 1, let f r = Jjf/ : [0,T] — > H be 
Lipschitz continuous and K, g W / f(IR + ) with \\K,\\ w[(R+) < 1- We also, recalling 
T>(A) = H 2 H V, assume the following compatibility conditions: 
forr = 2k(k = l,2,---), 

A k u° g V(A), A k u l g V, 

(2.16) f r -^(o)eH^nv, j = !,■■■ , k, k>l, 
^r-!W+i(o) g h*U-i) n V, j = 1, •• • , k, k > 1, 

and forr = 2k + l (k = 0, 1, 2, • • • ), 

A k u x g V(A), A k+1 u° g V, 

(2.17) / r " 2j '(0) Gff 2j ny, j = I,- - ,fe, k > 1, 

jr-2j+l^ g H 2U-1) n y 5 j = 1,. .. , * + 1, jfc > 0. 

T/ien £/iere is a unique solution of ( |1.1[ ). 

Moreover, for some C = C(k, \\IC\\w r (R+), T): 
for r = 2k (k = 1, 2, • • • ), we /ia/ye i/ie regularity estimate 

IK(t)|k + |K +1 (i)|| 

: r| ||A fe+1 u°|| + ||AV|| 



(2.18) 

fc ft „ t 



^pV''- 2j (o)||+^||^"- 1 / r - 2j+1 (o)||+ / ||/ r (s)||rfA 

,•-1 o-l JO / 



and, for r = 2fc + 1 (k = 0, 1, 2, • • • ), we /iave i/ie estimate 
IK(t)||v + ||u p+1 (*)ll 



(2.19) 



< C^||A fe+1 u || + || J 4 fe+1 u 1 || 

.7=1 .7 = 1 70 



Proof. 1. The case r = 1 follows from Theorem [T] Then, for a given r > 2, we 
show that 

(i) z r (0) = z r ' = (u r (0),u r+1 (0),0) g V(A), 

(ii) F r is differentiable almost everywhere on [0, T] and F g Li([0,T]; Z). 



CGFEM FOR HYPERBOLIC INTEGRO-DIFFERENTIAL EQS. 



9 



These imply existence of a unique strong solution of the abstract Cauchy problem 
(2.13), by [lis Corollary 4.2.10], that yields existence of a unique solusion of (1.1), 
by Lemma [2] 

2. First we note that (i) holds, if u r (0) G V(A) and u r+1 (0) G V, by the 
definition of T>(A). This can be verified by applying the compatibility conditions 
j2T6}-([2T7t in pl^- pl^ , using plo| . 

3. Now we prove (ii). By assumption f r : [0,T] — > H is Lipschitz continuous. 
Therefore, by a classical result from functional analysis, f r is differentiable almost 
everywhere on [0, T] and f r G Li([0,T]; H), since H is a Hilbert space. Then, 
recalling the assumption JC € W r f(M + ) and the fact that 

r 

f(i) = f(() + i^ u "(o)r-"((), 

n=0 

we conclude that F r (t) = (0, f r (t),0) is differentiable almost everywhere on [0,T] 
and F r G Li([0, T]; Z), that completes the proof of (ii). 

4. Hence, since A generates a Co-semigroup of contractions on Z by Corollary 
1, we conclude, by [TTJ Corollary 4.2.10], that there exists a unique strong solution 
z r = (u r ,v r 1 w r ) for the abstract Cauchy problem (2.13). This, by Lemma [2j 
proves that there is a unique solution u of (1.1), that completes the first part of 
the theorem. 

5. Finally, we prove the regularity estimates (2.18) and (2.19) for r > 2, since 
the case r = 1 follows from Theorem □ 

The unique strong solution z r — (u r ,v r ,w r ) of (2.13), is given by 

z r {t) = e tA z r ' n + f e^- s)A F r (s)ds, 



and we recall the fact that ||e*^||z < 1, since A is an infinitesimal generator of a 
Cq semigroup of contractions on Z. Therefore 



\\z r (t)\\z<\\z rfl \\z+ / \\F r ( S )\\ z ds. 
Jo 

Since v r = u r = u r+1 , z''<° = (u r (0), u r+1 (0), 0), and 

r-l 

\\F r (s)h - \\f r (s)\\ = \\f(s)\\ + £ \\Au n (0)\\\JC r - n -\ S )\ + \\Au r (0)\\\a S )l 

therefore we have 

((l- K )\\u'-(t)\\l + \\u r + 1 (t)\\ 2 ) j 

<((i- K )K(o)||^ + K +1 (o)|| 2 ) 1/2 



n=0 



1/2 



[]\r(s)\\ ds + 2ll^ n (0)ll /V - "- 1 ^)! ds 
Jo n=0 Jo 



\Au r 
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Hence, considering the assumption that ||JC||w^(r+) < 1, we have, for some C 
C(k, \\]C\\ W r T), 



(2.20) 



\\u r (t)\\ v + \\u r+1 (t)\\ < C[ \\u r (0)\\ v + K +1 (0)| 



\r( S )\\ d S +j2\\Au n (o)\\ 

n=0 



Since, by elliptic regularity estimate (2.6), 



so we have 



\\u r (0)\\ v <\\u r (0)\\ H ,<C\\Au r (0)\\, 

\\u r (t)\\ v + ik+^h < c(\\vr+Ho)\\ + E H^"(°)ll + /V( S )H 

\ „_n JO 



ds 



that, by ( |2.10[ )-( |2.12[ ), implies the regularity estimates ( |2.18[ )-( |2.19[ ). Now, the 
proof is complete. □ 

3. The spatial finite elment discretization 



The variational form of ( 1.1 ) is to find u(t) £ V, such that u(Q) — u°, u(0) = u 1 , 
and for t € (0,T), 

ft 



(3.1) 



(it, v) + a(u, v) 



IC(t — s)a(u(s), v) ds — (f, v), \fv £ V. 



Let n be a convex polygonal domain and {Th} be a regular family of triangula- 
tions of with corresponding family of finite element spaces Vfc C V, consisting of 
continuous piecewise polynomials of degree at most I — 1, that vanish on dfl (so the 
mesh is required to fit dQ). Here I > 2 is an integer number. We define piecewise 
constant mesh function Hk(x) — diam(if) for x £ K, K £ Th, and for our error 
analysis we denote h = m&XK^Th ^k- We note that the finite element spaces V h l 
have the property that 



(3.2) mm{\\v - X \\ +h\\v- X \\i} < Cti\\v\\ 



for v E H l n V, 1< i < I. 



We recall the ^-projection Vh '■ H — > Vl and the Ritz projection IZh ■ V — > V h l 
defined by 

a(TZ h v,x) = a(v,x) and {V h v,x) = {v,x)> V% € V^. 
We also recall the elliptic regularity estimate 



2.6 1, such that the error estimates 



(3.2) hold true for the Ritz projection TZh, see [2D]> i.e., 

(3.3) \\{K h -I)v\\ + h\\{K h - I)v\\! < Ch^vWu for v £ W n V, 1 < i < I. 



Then, the spatial finite element discretization of ( |3.1[ ) is to find Uh(t) £ such 
that u h (-, 0) = u° h , u h (-, 0) = u\, and for t £ (0, T), 

(3.4) {u h ,v h ) + a(u h ,v h ) - K(t- s)a(uh(s),Vh) ds = (f,v h ), Vv h £V, l l , 



where u° and u\ are suitable approximations to be chosen, respectively, for uq and 



u 1 in V^. 
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Theorem 3. Assume that 51 is a convex polygonal domain. Let u and be, 
respectively, the solutions of (3.1) and (3.4 1. Then 

(3.5) \\u h (T) - u(T)\\ <C\\ul-K h u Q \\ +Ch l (\\u(T)\\ l + \\u\\id?j, 

where we assume the initial condition u\ — VhU 1 . 

Proof. The proof is adapted from [4] . We split the error as 

(3.6) e = Uh — u = (v,h — Rhu) + {Rnu — u) = 9 + u>. 



We need to estimate 9, since the spatial projection error uj is estimated from (3.3). 
So, putting 9 in (3.4) we have, for vt € VJ[, 



, Vh) + a(9, Vh) - / K.(t - s)a(9(s),v h )ds 



(u h ,v h ) + a(u h ,v h ) - I K{t- s)a(uh(s),Vh)ds 



(R h u,v h ) - a(R h u,v h ) + / K(t— s)a(Rhu(s),Vh)ds, 



that, using (3.4), the definition of the Ritz projection TZh, and (3.1), we have 



v h ) + a(9, v h ) - I K,(t - s)a(9(s),v h )ds 



(3.7) 



(f,v h ) - (R h u,v h ) - a(u,v h ) + I K(t- s)a(u(s),v h )ds 
(u,v h ) - (R h u,v h ) = -(u>,v h ). 



Therefore we can write, for Vh(t) G V^, t e (0, T], 



dt 



,Vh) - {9,v h ) +a(9,v h ) - I IC{t- s)a(9(s),v h {t))ds = - — (uj,v h ) + (w,v h ), 



that, recalling e = 9 + uj, we obtain 



(3.8) - (9,v h ) +a(6,v h ) -J K(t - s)a(6{s),v h (t))ds = -^(e,v h ) + (u),v h ) 
Now let < e < T , and we make the particular choice 

V h (;t) = J 9{;T)dT, 0<t<T, 

then clearly we have 



(3.9) 



v h (-,e)=0, -V h (;t) = -9(;t), 0<t<T. 



Hence, considering (3.9) in (3.8) , we have 



Id 
2dt 



l^| 2 -!l^lly)- / lC(t-s)a(6(s),v h (t))ds = --(e,v h )-(uj,d). 
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Now, integrating from t = to t = e, we have 

r e r t 

2 \\a(n\\\2 ll„. /_M|2 i ||„, /nM|2 



||d( e )r-||fl(0)||' 1 -|K( E )||^ + |K(0)||^-2 / / K(t-s)a(6(s),v h (t))dsdt 

Jo Jo 

= -2(e(e) > « fc (e))+2(e(Q),u /l (0))-2 f (u,B)dt. 

Jo 

Then, using the initial assuption u\ = VhU 1 that implies the second term on the 
right side is zero and recalling Vh(e) = 0, we conclude 



re r t 

2 



(3.10) 



\WW + K(0)||f, -2 / / !C(t-s)a(9( S ),v h {t))dsdt 



JO 



< \\9(0)\\ 2 + 2 max ||0(t)|| / \\u>\\dt. 



Now, by changing the order of integrals, using 4r£(t — s) = —K{t — s) from (1.4 1, 
and integration by parts, we can write the third term on the left side as 



e rt re re 

(If 



-2 I I lC{t-s)a(6(s),v h {t))dsdt = 2 I I -f(t - s)a(6{s),v h (t))dtds 

JO JO JO Js 



2/ Z(e-a)a(0(8),v h {e))d8-2 J i{fS)a{9{a),v h {8))da 

2 So is ^ ~ S ^ a ( 9 ( S ")> i ' h ( t ^ dtdS - 



Then, using (3.9) and £(0) = K, we have 

re r t 

K{t ~ s)a{6(s) 1 v h {t))dsdt 



o Jo 



<\\v h {e)\\V -\\vh{Q)\\v) + 2 / / t{t-s)a{9{s),9{t))dtds 



Therefore, using this and Vh(e) = in (3.101 we have 

\\9(e)\\ 2 + (1 - K)||«fc(Q)||^ + 2 f a [ £(* - s)a(d(s),9(t))dtds 
< \\9(0)\\ 2 + 2 max \\9{t)\\ [ \\(lj\\dt, 

0<t<e J Q 

that considering the fact that £ is a positive type kernel and k < 1, in a standard 
way, implies that 

< c(||0(o)|| + CmW). 



Hence, recaling (3.6), we have 

||e(r)||<||0(T)|| + HT)|| 

< C(\\u° h - R h u°\\ + I \\u- R h ii\\dt) + \\(R h u - u)(T)\\, 



that using the error estimate (3.3) implies the a priori error estimate (3.5). □ 



CGFEM FOR HYPERBOLIC INTEGRO-DIFFERENTIAL EQS. 



13 



4. The continuous Galerkin method 
Here we formulate a continuous space-time Galerkin finite element method of 



order one, cG(l)cG(l), for the primar and dual problems (4.4 1 and (4.8 1, that is 
based on a similar method for the wave equation in [Jy. Then we prove stability 
estimaes for the discrete dual problem. These are then used in a priori error analysis, 
that is via duality. 



4.1. Weak formulation. First we write a "velocity-displacement" formulation of 



(1.1) which is obtained by introducing a new velocity variable. We use the new 
variables u\ — u, u 2 — u, and u — (iti,u 2 )> then the variational form is to find 
ux{t), u 2 (t) € V such that ui(0) = u°, u 2 (0) = v°, and for t € (0,T), 

(ui(t),«i) - (u 2 (t),ui) = 0, 

(4.1) (u 2 (t),v 2 ) +a(w 1 (t),u 2 ) - / lC(t - s)a(ui(s),v 2 ) ds 

Jo 

= (/(*), va), \f Vl ,v 2 eV. 

Now we define the bilinear and linear forms B : U x V — >• R and L : V — > R by 

B(u,v)= j - {u 2 ,v\) + (v,2,V2) + a(ui,v 2 ) 

JC(t — s)a(ui(s), v 2 ) ds > dt 



(4.2) 







(«i(0),vi(0)) + («2(0),« 2 (0)), 



T 



L{v)= / (/, V2 )di+(u , Ul (0)) + (« , W2 (0)), 



where 

U = H l (0,T;V) x ^(O.T; H), 
(4.3) r , / 

V = {v = («i,w 2 ) : to € L 2 (0,T;H) x L 2 (0, T; V), u.; right continuous in t). 



We note that the weak form (4.1 ) can be writen as: find u £M such that, 



(4.4) B(u,v)=L(v), V«eV. 

Here the definition of the velocity u 2 — iii is enforced in the L 2 sense, and the initial 
data are placed in the bilinear form in a weak sense. A variant is used in {7! where 
the velocity has been enforced in the H 1 sense, without placing the initial data in 
the bilinear form. We also note that the initial data are retained by the choice of 
the function space V, that consists of right continuous functions with respect to 
time. 

Our a priori error analysis for the full discrete problem, cG(l)cG(l) method in §6, 



is based on the duality arguments, and therefore we formulate the dual form of (4.4 ) 



To this end, we define the bilinear and linear forms B* : V* xW* — > R, L* : V* 
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for r € by 

B*(v, z) — J | — (ui,ii) +a(vi,Z2) — J fC(s - t)a(vi,z 2 (s)) ds 
(4.5) - {v 2 ,z 2 ) - (v 2 , Zl )}dt+ («i(T),zi(T)) + (u 2 (T),z 2 (r)), 

L»= / {( Ul! j 1 ) + ( U2! j 2 )}dt+( Ul (T),zn + (^(T),z 2 T ), 



where j 1; j 2 and z^, zj represent, respectively, the load terms and the initial data 
of the dual (adjoint) problem. In case of r = 0, we use the notation B*,L* for 
short. Here 

U* = H\0,T;H) x H l (Q,T;V), 
(4.6) r 

V* = {v = (vi,v 2 ) G L 2 (0,T;V) X L 2 (0,T; H) : v t left continuous in t}. 



We note that, recalling tyL3$, U C V*, U* C V. 

We also note that B* is the adjoint form of B. Indeed, integrating by parts with 
respect to time in B, then changing the order of integrals in the convolution term 
as well as changing the role of the variables s, i, we have, 

(4.7) B(u,v) = B*(u,v), VueU,veU*. 

Hence, the variational form of the dual problem is to find z *E U* such that, 

(4.8) B*{v,z)=L*{v), VweV*. 

4.2. The cG(l)cG(l) method. Let = t < t x < ■ ■ ■ < t n ^ 1 <*„<•••< 
t-N = T be a partition of the time interval [0, T]. To each discrete time level t n we 
associate a triangulation 7^™ of the polygonal domain f2 with the mesh function, 

(4.9) h n {x) = h K = diam(i-Q, x € K, K e 7?, 

and a finite element space VP consisting of continuous piecewise linear polynomials. 
For each time subinterval /„ = (t n —i,t n ) of length k n — t n — t n -i, we define 
intermediate triangulaion 7^ ra which is composed of the union of the neighboring 
meshes 7^", 7^ l_1 defined at discrete time levels t„, f„_i, respectively. The mesh 
function h n is then defined by 

(4.10) h n {x) =h K = diam(K), x G K, K G f£. 

Correspondingly, we define the finite element spaces Vfi consisting of continuous 
piecewise linear polynomials. This construction is used in order to allow continuity 
in time of the trial functions when the meshes change with time. Hence we obtain a 
decomposition of each time slab f2 n = ft x I n into space-time cells K x I n , K £ 7^" 
(prisms, for example, in case of ft C M 2 ). The trial and test function spaces for the 
discrete form are, respectively: 

Uhk = \u = {U\,U 2 ) : U continuous in Vt x [0,T], U(x,t)\i n linear in t, 

(4.11) J 

V hk = {V= (Vi,V 2 ) : V{-,t) continuous in fl,V{-,t)\ In G (V^ 1 ) 2 , 



V(x,t)\i n piecewise constant in t|. 
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We note that global continuity of the trail functions in Uhk requires the use of 
'hanging nodes' 1 if the spatial mesh changes across a time level t n . We allow one 
hanging node per edge or face. 

Remark 3. If we do not change the spatial mesh or just refine the spatial mesh 
from one time level to the next one, i.e., 

(4.12) V^CV,:, n=l,...,N, 
then we have V£ = V£. 

In the construction of Uhk and Vhk we have associated the triangulation T h n with 
discrete time levels instead of the time slabs Vl n , and in the interior of time slabs we 
let U be from the union of the finite element spaces defined on the triangulations 
at the two adjacent time levels. This construction is necessary to allow for trial 
functions that are continuous also at the discrete time leveles even if grids change 
between time steps. For more details and computational aspects, including hanging 
nodes, see |14j and the references therein. Associating triangulation with time slabs 
instead of time levels would yield a variant scheme which includes jump terms due 
to discontinuity at discrete time leveles, when coarsening happens. This means 
that there are extra degrees of freedom that one might use suitable projections for 
transfering solution at the time levels t n , see [7]. 

The continuous Galerkin method, based on the variational formulation (4.1), is 
to find U € Uhk such that, 

(4.13) B(U,V)=L(V), We V hk - 
Here, as a natural choice, we consider the initial conditions 

(4.14) ul := E7i(0) = V h u°, u\ := U 2 (0) = V h u\ 



where the L 2 pojection Vh is defined in (4.201. 

The Galerkin orthogonality, with u = (1*1,1(2) being the exact solution of (4.1), 
is then, 

(4.15) B(U-u,V)=0, W eV hk - 

Similarly the continuous Galerkin method, based on the dual variational formula- 
tion (4.8), is to find Z 6 Uhk such that, 

(4.16) B*(V,Z)=L*(V), We V hk , 
Then, Z also satisfies, for n = 0, 1, . . . , N — 1, 

(4.17) BUV,Z) = LUV), W eVhk- 



From (4.13) we can recover the time stepping scheme, 
f {(U 1 ,V 1 )-(U 2 ,V 1 )}dt = 0, 

f {{U2,V 2 ) + a{U 1 ,V 2 )- J^K{t-s)a{U 1 {s),V 2 )ds]dt 
= [ {f,V 2 )dt, V(V U V 2 ) e Vhk, 



(4.18) 



U 1 (0)=<, U 2 {Q)=ul 
with the initial conditions ( 4.14[ ). 
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Typical functions U = {U\, U 2 ) G V hk , W = (W 1 ,W 2 ) G Whk are as follows: 

U i (x,t n ) = U?(x) = J2ur ij ^(x), 

3=1 

(4.19) Ui{x,t)\ In = ^ n . 1 {t)U^- 1 {x)+i) n {t)U^{x), 

Wi(x,t)\ In =J2W^(x), 

3=1 

where m n is the number of degrees of freedom in {f" (x)}J^i are the nodal 
basis functions for V£ defined on triangulation 7^™, and ip n {t) is the nodal basis 
function defined at time level t n . Hence (4.18) yields 

M njjn _ _^ M njjn = ^-n-l.n^n-l + > 

M n u% + - w- n )s»u? = A-r-^ur 1 + (- y + c+js-^c/r 1 

n-1 



where, for I = 1, . . . , n 



tAti r pt/\ti 

u ti.i = I I JC(t — s)i/)i-i(s) ds dt, = / / JC(t — s)i/>t(s) ds dt, 



s ,,n = (s!f)ij = K^,v?)) iJ . 

We define the orthogonal projections 7£/i, n : — > V^", 7\ n : H V£ and 
P*,n : ^ 2 (4) d -> Po(^„), respectively, by 

a(K Kn v -v, X ) = 0, Vv G V, x G V?, 
(4 20) (V hin v-v, X ) = 0, VveH, x eV£, 

J{V k .nV -v)-i/)dt = 0, toe L 2 (I n ) d , i> e Pg(I n ) , 

with Pq denoting the set of all vector-valued constant polynomials. Correspond- 
ingly, we define TZhV,VhV and Vkv for t € I n (n = 1> - "" ,N), by (JZhv)(t) — 
K hjn v(t), (V h v)(t) = Vh, n v(t), and V k v = 7\„(u|/J. 

Remark 4. In the case of assumption (4.12), by Remark [3] and the definition of 
the L 2 -P r ojection Vk, we have V, V^V 6 Vhk, for any V G Uhk- 

We introduce the discrete linear operator A n:T : V£ — > V£ by 
a(v r ,w n ) = (A nir v r ,w n ), to r e V£, w n e V£. 
We set A„ = A„ „, with discrete norms 



\\v n \\ h ,l = ll^ /2 «nll = \j(v n ,A l n v n ), v n e V£ and / G 
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and Ah so that Ah,v — A n v for v € V,™. We use Ah when it acts on Vfi. For later 
use in our error analysis we note that V%A — AhTZh- 



4.3. Stability of the solution of the discrete dual problem. We know that 
stability estimates and the corresponding analysis for dual problem is similar to the 
primar problem, however with opposite time direction. Hence, having a smooth or 



weakly singular kernel with (1.2), we can quote slightly different energy identities, 



compare to (4.21), from [7] or [3] for the discrete dual solution, from which simi- 



lar stability estimates to (4.22) is obtained, though with different projections and 
constants. 

To prove stability estimates in [7] and [16] we have used auxiliary functions in 
the form, respectively, W(t, s) = U{t) - U{s) and W(t, s) = U(t) - U(t - s). Here, 
using the properties of the fuction £ = in the convolution integral and partial 
integration, we give a proof which is straightforward. 



We note that the stability constant in (4.22) does not depend on t. See [13] . [T5] 
and [12] , where stability estimates have been represented, in which the stability 
factor depends on t, due to Gronwall's lemma. 



Theorem 4. Let Z be the solution of (4.16) with sufficiently smooth data zf,Z2, 
ji,j2- Further, we assume (4.12). Then for ! e R, we have the identity, 



+ k\\Z 2 (t n )\\l l+1 +2^ J £(s - t)a(A l h Z 2 (t),Z 2 (s)) dsdt 
= \\Z 1 (T)\\l d + {l + K)\\Z 2 {T)\\l, +1 



(4.21) 



2^ {A l h Z l ,V k V h j 1 )dt + 2 jf (A l h +1 Z 2 ,V k r h j 2 )dt 

2 J J IC{s~t)a(A l h V k V h j 2 ,Z 2 { s ))dsdt 

-2 IC(T-t)a(A l h Z 2 (t),Z 2 (T))dt 
-2£(T-t n )a(A l h Z 2 (t n ),Z 2 (T)), 
where h= 1 — k. Moreover, for some constant C — C(k), we have stability estimate 
\\Zi(t n )\\h.i + \\Z 2 (t n )\\ h>l+1 < C^Z^T)^ + ||Z 2 (T)|| htJ+1 

V h ji\\ h ,i + \\V h j 2 \\ h ,i+i) dt}. 



(4.22) 



Here, we set the initial data of the dual problem as 



(4.23) 



Zi(T) = V h zJ, i = 0,1. 



Proof. 1. The solution Z of ( |4.16| ) also sati sfies ( |4.17[ ), for n = N - 1, 
Then recalling Remark [4] for the assumption (4.12), we obviously have, 



,1,0. 



(4.24) 



v k Zy = -z 2 - r k r h j 2 . 
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2. Using this in (4.171 we obtain 

cT 



(V X ,Z X ) + a(V x , Z 2 ) -J K{s- t)a(V x ,Z 2 (s)) ds} dt + (Vi(T), Z X (T)) 
+ (V 2 (T),Z 2 (T)) = f (V x ,j x )dt + {V x {T),zJ) + (U 2 (T),z 2 T ). 



For the convolution term we recall K.(s — t) = — Z? s £(s — t) from ( |1.4[ ), and then 
partial integration yields 



T P T 



T r T 



/C(s - t)a(Vi, Z 2 {s)) ds dt 

£(s - t)a(V x ,Z 2 {s)) dsdt + jf £(T - t)a(Vi, Z 2 (T)) dt 

-k£ a(Vi,Z 2 (t))dt. 
These and k = 1 — K imply that the solution Z satisfies, 

(V x ,Z x ) + ka(V x ,Z 2 )- J Z(s-t)a(V x ,Z 2 (s))ds 

+ £(T - t)a(V x ,Z 2 (T))} dt + (Vi(T), Z X {T)) + (V 2 (T), Z 2 {T)) 

T 

{Yu'Phjx) dt + (V X {T), zj) + (V 2 (T),z^), 



Now we set Vi = A l h V k Zi, and recall the initial data (4.231 such that the terms con- 
cerning the initial data are canceled by the definition of the orthogonal projection 



Vh- Then we have 
(4.25) 

T , 

ill 



-{A l h V k Z x ,Z x ) + m(A l h V k Z x ,Z 2 )- J £(s - t)a(A l h V k Z x , Z 2 {s)) ds 

+ f (T - t)a(A l h V k Z x , Z 2 (T))} dt = jf (A l h V k Z x ,V h j x ) dt. 

3. We study the four terms at the left side of the above equation. For the first 
term we have 



(4.26) 



-(A l h V k Z x ,Z x ) dt = -±\\Z x {T)f h!l + \\\Z x (t n )\\l,- 



With (4.24) we can write the second term as 

i-T 



is J a{A l h V k Z u Z 2 )dt 



(4.27) 



= — « 



a(A l h Z 2 , Z 2 ) dt-R a(A l h V k V h j 2 ,Z 2 ) dt 



;\\MT)\\il+i + ^\\Z2(t n )\\li +x -K a(A l h V k V h32 ,Z 2 )dt 
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For the third term in (4.25), by virtue of ( 4.24 1 and integration by parts, we obtain 

cT r T 



£(s - t)a(A l h P k Z 1: Z 2 {s)) dsdt 



(4.28) 



T fT 



£(s - t)a(A l h Z 2 {t), Z 2 (s)) ds dt 



T r T 



t„ Jt 
T 



K{s - t)a(A l h V k V h j 2 , Z 2 (s)) dsdt 



£(T - t)a(A l h V k V h j 2 , Z 2 {T)) dt- k a(A l h V k V h j 2 , Z 2 {t)) dt 



Finally, for the last term at the left side of (4.251, we use (4.24) and integration by 
parts to have 



J £(T-t)a(A l h V k Z 1 ,Z 2 (T))dt 



(4.29) 



K(T - t)a(A l h Z 2 (t),Z 2 (T)) dt - «||Z 2 (r)||L +1 



+ £(T - t n )a(A l h Z 2 (t n ), Z 2 (T)) £(T - t)a(A l h V k V h j 2 , Z 2 {T)) dt. 



Putting ( |4.26[ )-( |4729] ) in ( |4.25| ) we conc lude the identity ( |4.21 |. 

4. Now we prove the estimate (4.22). We recall, from (1.5), that £ is a positive 
type kernel. Then, using the Cauchy-Schwarz inequality in (4.21 ) and ||£||li(r+) — 
k, £(i) < k, we get, for C 3 = C 3 (k) and C 4 = C^k), 

\\Mtn)\\h + m 2 (t n )\\i l+1 

< ||^(T)||^ + (i + «)||z 2 (r)||^ +1 

+ d ( max T HZill^ + 1/Ci( jf \\P h V h h\\ hl l dt)" 

+ C 2t razxJZ 2 \\l l+1 + l/C 2 (J \\P k P h j 2 \\ h ,i + i dt) 2 
+ C 3 1 1 Z 2 (T) + 1/C 3 t max, 1 1 Z 2 \ \ 2 hd+1 
+ C4Z 2 (T)\\l l+1 + l/C4Z 2 (t n )\\l l+1 . 
Using that, for piecewise linear functions, we have 



(4.30) 
and 



max \Ui\ < max \Ui(t n )\, 

[0,T] 0<n<N 



/ \V k f\dt< f \f\dt, 
Jo Jo 

and that the above inequality holds for arbitrary TV, in a standard way, we conclude 



the estimate inequality (4.22). Now the proof is complete. 



□ 
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5. A PRIORI ERROR ESTIMATION 



We define the standard interpolant I k with I k v belong to the space of continuous 
piecewise linear polynomials, and 

(5.1) I k v(t n ) = v(t n ), rc = 0, 1, ••• , N. 

By standard arguments in approximation theory we see that, for q = 0, 1, 



\\I k v-v\\idt < Ck q+1 / \\D% +1 v\\idt t for i = 0, 1, 



(5.2) 



where k = maxi<„<jv k n . 

We recall that we must specialize to the pure Dirichlet boundary condition and a 
convex polygonal domain to have the elliptic regularity (2.6 1, from which the error 
estimates (3.3) hold true for the Ritz projections in (4.20|. We note that the energy 
norm || • \\y is equivalent to || • |ji on V. 



Lemma 3. Assume (4.12). Then, for V, W £ Uy lk , we have 



(5.3) 



B*(V k V,W) = B(V,V k W) 

+ (Vi(0), (Wi - 7Wi)(Q)) + (V 2 (0), (W 2 - 7W 2 )(Q)) 

- (Vi(T), Wi(T)) - (V 2 (T),W 2 (T)) 

+ ({V k Vi)(T),W!(T)) + ((P k V 2 )(T),W 2 (T)). 



Proof. We recall Remark 4 for the assumption ( |4.12[ ), and the definition of the 
bilinear forms B,B* from (4.2) and (4.5). Then by the definition of V k and partial 
integration in time we have 

B*(V k V,W) 



- (V 2 ,W 2 ) - O^TWr)} 



K,{a-t)a{V 1 ,V k W 2 {s)) ds 



dt 



((Wi)(T),WMT)) + ((V k V 2 )(T),W 2 (T)) 

{(Vl,Wi) + a(V u P k W 2 ) - J K(s-t)a{VuV k W 2 (s)) ds 

(V 2 ,W 2 ) - (V 2 ,7Wi)} dt 
+ (^(0),W!(0)) + (V 2 (0),W 2 (0)) - (Vi(T),Wi(T)) - (V 2 (T),W 2 (T)) 
+ {(VkVt^WiiT)) + ((P k V 2 )(T),W 2 (T)), 
that implies (|5.3|), and the proof is complete. □ 



Theorem 5. Assume that f2 is a convex polygonal domain, and ( 4.12[ ). Let u and 

U be the solutions of (4.4) and (4.13). Then, with e = U — u and C = C(k), we 
have 



(5.4) 



ei(T)||<C/> 2 || U || 2 + K(T)|| 



"ill 2 dt 



Ct 



(lltiall + llfiiHi)*, 
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and, with a quasi-uniform family of triangulations, 

\\e 1 (T)\\ 1 <Ch(\\u°\\ 1 + \\u 1 (T)\\ 2 + f \\inhdt) 
(5.5) V 



+ Ck 2 (Hfialli + lluilla)*, 



\\e 2 (T)\\<Ch(\\u% + \\u 2 (T)\\ 



\ui\\ 2 dt 



+ Ck 2 (\\U2\U + \\ Ul \\ 2 ) dt. 



(5.6) 



Proof. 1. We recall Remark [4] for the assumption ( |4.12[ ). We set 

(5.7) e = U — u = (U — h^hu) + (Ik^hU — n h u) + (n h ,u — u) = 9 + 77 + u, 

where I k is the linear interpolant defined by ( |5.1[ ), and 717, is in terms of the pro- 
jectors IZh and Vh, such that 

6»i = Ui - I k TZ h ui, rji = (h - I)7l h ui, oji = (U h - I)ui, 



(5.8) 



h = U 2 - IkVhU 2 , m = (h - I)Vh.u 2 , lo 2 = (V h - I)u 2 . 



We note that 77 and u> can be estimated by ( |5.2[ ) and (3.3), and therefore we need 
to estimate 9. 

2. Now, putting V — V k 9 in (4.16) with j\ = j 2 = 0, we have 
(5.9) L*(V k 9) = ((M)(n^ T ) + ((^ 2 )(T),z 2 T ) = B*(V k 9,Z), 
that, using Lemma [3] and the initial data ( |4.23[ ), implies 
(9 1 (T),Z 1 (T)) + (9 2 (T),Z 2 (T)) 

= B(9,V k Z) + (^(0), (Z 1 - V k Zi)(0)) + (&(Q), (Z 2 - V k Z 2 )(0)) 



{OuVkZx) - (0 2 ,'P fc Zi) + (6 2 ,V k Z 2 ) + a{6i,V k Z2) 



K(t - s)a(9 1 (s),V k Z 2 ) ds\ dt 
(ei(O).Zi(O)) + (0 2 (O),£ 2 (O)). 



Then, using 9 = e — 77 — lo and the Galerkin orthogonality (4.15), we have 
(9 1 (T),Z 1 (T)) + (9 2 (T),Z 2 (T)) 

= J Q {- faVkZi) + {r)2,V k Zi) - (r) 2 ,V k Z 2 ) - a{ m ,V k Z 2 ) 

+ J fC(t- s)a( m {s),V k Z 2 )ds}dt 
-( m (0),Z 1 (0))- (772(0)^2(0)) 
+ J o {- fatPhZi) + {oj 2 ,V k Zi) - {d) 2 ,V k Z 2 ) - a(cj u V k Z 2 ) 

+ J fC(t- s)a(wi(s),V k Z 2 )ds}dt 
-(wi(0),Zi(0))-(wa(0),Z a (0)). 
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By the definition of rj, that indicates the temporal interpolation error, terms in- 
cluding rji and 7^(0) vanish. We also use the definition of u; in (5.8), that indicates 
the spatial projection error, and we conclude 

(e 1 (T) 7 Z 1 (T)) + (9 2 (T) 7 Z 2 (T)) 

fatVkZi) - a{r) U V k Z 2 ) + I K(t - s)a(r ]1 (s),r k Z 2 ) ds} dt 



{Cj 1 ,V k Z 1 )dt-(uj 1 (0),Z 1 (0)), 



that, setting the initial data Zx{T) = V h z[ = A^ l 6i(T) and Z 2 {T) = V h z 2 = 



A] 



' 2 (T), / G K and using the Cauchy-Schwarz inequality, we have 
Pi(T)\\l_ l + \\9 2 (T)\\l_ (l+1) 

<dnax\\PkZl\\h + Vd( [ WmWK-idtf 

— — J 



(5.10) 



+ d max \\V k Z 2 \\ 2 hl+1 + l/d 

0<t<T !+ 



+ C 3 max \\V k Z 2 \\ 2 hl+1 + l/d 

0<t<T ii/i.i + i 



+ d max \\V k Zi\\ 2 hl + l/d 

0<t<T ' 



T > 2 
||?7l||fc_; + l dt . 



(lC*\\Th\\h,-l+l){t)dt 



\\V h u x \\h,-idt\ 



+ C B ||Zi(0)||^ + l/C4V h ui (0)\\l-i- 

On the other hand, putting the initial data Z X (T) = A^ l 6i(T) and Z 2 (T) 
A h {l+1) 6 2 {T) in the stability estimate ( |4.22[ ) with j x = j 2 = 0, we have 

||Si(*»)IU,i + \\z 2 {t n )\\hj+i < c{||»i(r)lU,_i + \\e 2 (T)\\ ht _ (l+1) }. 

Using this, together with (4.30), and ||£||l 1 (k+) 
we have 



k in (5.10), in a standard way, 



pi 

< c{\\V h wi(0)\\h,-i + / (lMk-i + IMk-i+i + HPhWilk-j) dt}. 



(5.11) 



3. To prove the first error estimate (5.4), we set I = 0, and we recall the facts 
that || • || M = || • ||, || • || M = || • ||i. Then, recalling e(T) = 9{T) + V (T) + lu(T) = 
0(T) + lu(T), and L 2 -stability of the projection Vh, we have 

\\ei(T)\\ < c[\\{K h - I)u°\\ + \\(K h - I)ux(T)\\ 

+ f (||(/fc-/)ti 2 || + ||(/fc-i)ti l ||i + ||(^/,-i)6i||)dt}. 



This completes the proof of the first a priori error estimate ( |5.4[ ) by (5.2) and (3.3). 

4. Now, to prove the last two error estimates (5.5)-(5.6), we set I = —1 in (5.11 1, 
and we recall the assumption of having a quasi-uniform family of triangulations. 
Then if -^stability of the I^-projection Vh, that is 



(5.12) 



im«||i<C|H|i, veH\ 
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holds true. Hence, recalling e(T) = 0(T) +u>(T), we have 
||e x (r)IU < C{\\(1l h - JKlli + - /)«i(T)|| x 

+ J (||(4 - /)«a||i + || (A - /)«i|| a + \\(K h - I)ui\\t) dt} 
\\e 2 (T)\\ < c{\\(K h - I)u% + \\(K h - I)u 2 (T)\\ 

+ [ (||(/j t -J)« 2 ||i + ||(Jfc-/)ui||a + ||(^fc--0«i|| 1 )dt}. 



This completes the proof of the error estimates (5.5)-(5.6) by (5.2) and (3.3). Now 
the proof is complete. □ 



We note that the assumption of quasi-uniformity for validity of (5.12), that is 
used for error estimates (5.5)-(5.6), can be relaxed, see [B|, though it is not an 



considerable restriction in a priori error analysis. 

6. Numerical example 

Here we verify the order of convergence of the cG(l)cG(l) method by a simple 
example for a one dimensional problem with smooth convolution kernel. Another 
example for two dimensional case with similar results, with fractional order kernels 
of Mittag-Lefner type, can be found in [7]. 

We consider a decaying exponential kernel with ||/C||l 1 (k+) = k = 0.5, the initial 
data u° — u 1 = 0, and load term / = 0. We set homogeneous Dirichlet boundary 
condition at x = and a constant Neumann boundary condition at the end point 
x = 1, toward negative y axis. Figure 1 shows that the method preserves the 
behaviour of the model problem. 

In Figure 2, we have verified numerically the spatial rate of convergence 0(h 2 ) for 
L 2 - n orm of the displacement. In the lack of an explicit solution we compare with a 
numerical solution with fine mesh sizes h, k. Here h m i n — 0.0078 and k m i n — 0.017. 
The result for temporal order of convergence, 0(k 2 ), is similar. 



U(0.25,t) 

U(0.5,t) 

U(1 ,tj 




Figure 1. Damping of the oscillating material: at points x = 0.25, 0.5, 1. 
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Figure 2. Order of convergence for spatial discretization 
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